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^ ■ Abstract 

Wc introduce a nonsymmctric real matrix which contains aU the information that 
the usual Hermitian density matrix does, and which has exactly the same tensor 
^ ' product structure. The properties of this matrix are analyzed in detail in the case 

' of multi-qubit (e.g. spin = 1/2) systems, where the transformation between the 

■ real and Hermitian density matrices is given explicitly as an operator sum, and 

cn 
o 



X 



used to convert the essential equations of the density matrix formalism into the 



' real domain. 

^3 , PACS numbers: 03.65.Ca, 03.67.-a, 33.25.+k, 02.10.Xm 



I. PROLOGUE 

The density matrix plays a central role in the modern theory of quantum mechanics, and 
an equally important role in its applications to optics, spectroscopy, and condensed matter 
physics. Viewed abstractly, it is a self-adjoint operator p on system's Hilbert space, the ex- 
pectation values 0<(?/'|p|^/;)<lof which give the probability of observing the system in 
the state As a matrix, however, it is generally represented versus the operator (or "Li- 
ouville") basis induced by a choice of a complete orthonormal basis {|z)|z = 0,l,...} 

in the underlying Hilbert space. These complex-valued matrices p = [{i \p\j are neces- 
sarily Hermitian and positive semi-definite. Their diagonal entries are the probabilities of 
these mutually exclusive basis states, whereas their off-diagonal entries prescribe the amounts 
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by which the probabihties of their coherent superpositions deviate from the corresponding 
classical mixtures due to interference. 

Another option is to use a operator basis the elements of which have rank exceeding one, 
so that it is not induced by any Hilbert space basis. The most common example here is the 
representation of operators on a two-dimensional Hilbert space by real linear combinations 
of Pauli matrices {cro(= 120)7 <^i: ^s}- In this case the basis elements themselves are self- 
adjoint and so can be given a physical interpretation, e.g. as the components of the Bloch or 
Stokes vector (Bloch, 1946; Feynman et al., 1957). Although arbitrary bases of self-adjoint 
operators could be used, for multi-particle systems it is desirable that the overall basis be 
induced by identical bases on each particle's Hilbert subspace. With the Pauli matrices this 
leads to the so-called product operator representation (Ernst et al., 1987), in which density 
operators are represented by linear combinations of all possible tensor [ergo Kronecker) 
products of the Pauli matrices, e.g. crlcrj = {crj.^ (Tq){(Tq® cr^). In contrast to the Hermitian 
case, it has not been widely recognized that this tensor product structure is reflected by 
the real-valued coefficients in the expansion of any density operator in terms of product 
operators (also known as the coherence vector (Mahler and Weberruss, 1998)). Thus if one 
properly arranges these coefficients in a matrix one obtains a real but nonsymmetric analog 
of the Hermitian density matrix with the same tensor product structure. This fact holds for 
any number of Hilbert spaces S)k (/c = 1, 2, . . .) of arbitrary (even infinite) dimension L > 
and self-adjoint bases {-B^ j^Li for the space of bounded linear operators on each. 

The purpose of this paper is show how, in the case of multi-qubit {ergo two-state quantum) 
systems, one can perform all the usual Hermitian density matrix calculations entirely with 
the real density matrix. While the formulae are more complicated in most cases than they 
are with the Hermitian density matrix, we argue that they are in many respects closer 
to the underlying physics than the Hermitian formulae, simply because the entries of the 
real density matrix correspond to (expectation values of) observables. Indeed, it is well- 
known that the single-qubit Pauli algebra is nothing but a complex matrix representation 
of the geometric (or Clifford) algebra of a three-dimensional Euclidean vector space (Havel 
and Doran, 2002a). This real algebra in turn has been demonstrated to be a concise but 
versatile formalism within which to analyze and teach much of modern physics (Doran and 
Lasenby, 2003; Hestenes, 2003). It makes a certain amount of sense to use a representation 
wherein the "reverse-even" entities in the algebra, i.e. scalars and vectors, are real while the 
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less-familiar "reverse-odd" entities, i.e. "pseudo-scalars" and "pseudo- vectors" , are purely 
imaginary (cf. Baylis (1999)). The drawback of our representation is that the matrices no 
longer form a representation of the underlying geometric algebra, i.e. the geometric product 
no longer corresponds to matrix multiplication. We will leave it to the community to decide 
if or when the advantages outweigh the disadvantages, and offer our results simply as the 
outcome of an intellectual exercise. 

II. METAMORPHOSIS 

We begin by introducing a bit of notation which will considerably simplify the remainder 
of our presentation. First, instead of the above bra-ket notation, let us write the 2x2 
elementary matrices as 

Eoo^|0)(0|, Eio^|l)(0|, Eoi^|0)(l|, En^|l)(l|, (1) 

where |0) ^ eg, |1) ^ ei denote an orthonormal basis for a two-dimensional Hilbert 
space. Then it is easily seen that, for any nonnegative integers z,j < M = 2^ — 1, the 
(M -|- 1) X (M -|- 1) elementary matrix Ejj is the Kronecker product of 2 x 2 elementary 
matrices the indices of which are the bits in,jn ^ {0, 1} in the binary expansions of i,j, 
respectively, i.e. 

Eij = [SikSji] = Ei^j^ ® ■ ■ ■ Ei^j^ , (2) 

where the 5's are Kronecker deltas. In an analogous fashion, we will denote the usual 2x2 
Pauli matrices by 

Poo = CTo , Pio = CTi , Poi = (72 , Pn = (Ta . (3) 

In this notation it may readily be verified that the multiplication table among the Pauli 
matrices may be expressed succinctly as 

P,,P,, = ^(^^-^■^)(^-2^^')(^-2^^)P(,+,_2.fc),o+,-2i£^ (z,j,A;,£g{0, 1}) , (4) 

where = —1. This indexing scheme may be extended to all Kronecker products of these 
matrices in the same way as for the 2x2 elementary matrices, i.e. 

Pi, = P,,,, ® ■ ■ ■ ® P,^,^ (0 < 2, J < M) . (5) 
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For example, if M = 3 we have Pqi = ctq ® cr2, P02 = ® ctq, and P03 = (T2® (72- 

The Hermitian density matrix, of course, can be expanded relative to either the elemen- 
tary matrix basis or the Pauli matrix basis, e.g. 

1,1 1 ^-^1,1 



(6) 



for a single qubit with pij G C but cr^ 



Both of these bases are orthogonal relative 



to the Hilbert-Schmidt inner product, which is given by (X| Y) = tr(XY^) = 2^ (XY^) 
for any number of qubits > 0. Thus there is a unique unitary superoperator 
which carries the Pauli to the elementary matrix basis, where the factor of V2 comes from 
||Pjj|p = {Pij I Pij) = 2 ||Ejj|p. This superoperator, moreover, is nearly self-adjoint since 
( Eij I Pke ) = { Efcf I Pij ) for all < i,j,k,i < 1 with the sole exception of ( Eio | Pqi ) = 
z { Eqi I Pio ) and its complex conjugate. On applying U to both sides of Eq. (6), therefore. 



we obtain 
Uip) 



■1) 



fx- • E- • 



(7) 



^i,j=0 V ^ ' -i— /jj'=0 

We will take this as our definition of the real density matrix for a single qubit. Henceforth, 
we shall denote this by 



o"io cm 



Uip) 



(pCTo) (pO-2) 



(8) 



Note that our choice of normalization gives aoo = 2( p) = 1, so that although U is otherwise 
unitary the Hilbert-Schmidt norm is scaled by a factor of \/2; specifically 2 ||p|p = ||cr|p = 



01 



(9) 



Let us now look at some explicit representations of the superoperator U. To begin with 
the mapping from the Pauli to the basis elementary is clearly 

Eoo = KPoo + Pii) Eoi = i(Pio + ?Poi) 
Eio = |(Pio — «Poi) Ell = |(Poo ~ Pii) 
and hence (since coordinates are contravariant) 
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VW\p) ^ U\p), (10) 
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where we have factored the overall superoperator's matrix tl into the product of a diagonal 
matrix V and a purely real one W. An operator sum representation for the superoperator 
W may be derived from the singular value decomposition of its Choi matrix, i.e. 
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where the Choi matrix (left) is obtained simply by swapping certain pairs of the entries in 
W = [^ijjfjlo (Havel, 2003); specifically W20 ^ u^oi, ^730 ^ u^n, W22 ^ wqs and W32 ^ w^is- 
Observe that the left-singular vectors associated with the nonzero singular values of W can 
be written as "columnized" Pauli matrices, specifically | Pqo ) and | Pio ) in the notation of 
Eq. (10), while the corresponding right-singular vectors are | Pn ) and | Pio )• It follows 
that the operator sum form of W is (Havel, 2003, Proposition 3) 



W(p) = PoopPii + PiopP 
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Poo + pu Pio - Poi 
Pio + Poi Poo - Pll 



(12) 



Although an operator sum form for V could be obtained by this same approach, since it is 
diagonal a more compact representation of its action may be obtained by packing its nonzero 
entries into a single 2x2 matrix which acts via the "entrywise" or Hadamard product "0" 
(Havel et al., 2001), as follows: 



CT = Q© (pPii + PiopPio) = 



1 -I 
1 1 







Poo + Pll Pio - Poi 
Pio + Poi Poo - Pll 

Poo + Pll ^ (Poi - Pio) 
Poi + Pio Poo - Pll 



(13) 



Since W is self-adjoint and the overall superoperator hi is unitary (up to a factor of -\/2), it 
is easily seen that the inverse can be written as 



p = i((Q0o-)Pn + Pio (Q Oct) Pio) , 
where the overbar indicates the complex conjugate of all the matrix entries. 



(14) 
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The beauty of this operator sum form for the superoperator lA is that the Hadamard 
product obeys the mixed product formula with the Kronecker product, 

(A ® B) (C ® D) = (A C) ® (B © D) , (15) 

just hke the usual matrix product does. Thus if we extend lA to factorable multi-qubit 
density matrices in the obvious way, 

U{p^ ®---®p^) = U{p^) ® ■ ■ ■ ® W(p^) , (16) 

and thence to arbitrary multi-qubit density matrices by linearity, we immediately obtain a 
general expression. Explicitly, in the case of two qubits, we get 

= (qq {p" Pn + Pio Pio)) ® (q (p' Pii + Pio Pio)) 
= (Q Q) ((p' Pn + Pio p' Pio) ® (p' Pii + Pio p' Pio)) 

= (Q Q) ((p' Pn) (p' Pn) + (p' Pn) ® (Pio p' Pio) ■ ■ ■ ^^^^ 

+ (Pio p" Pio) ® (p' Pn) + (Pio Pio) ® (Pio Pio)) 

= (Q Q) ((p^ p')(Pii Pn) + (Poo ® Pio)(p' ® p')(Pn ® Pio) ■ ■ ■ 

+ (Pio Poo)(p' ® p')(Pio ® Pn) + (Pio ® Pio)(p' ® p')(Pio ® Pio)) 
= p^^ ^ p^^ p p^^ ^ p^^ p p^^ ^ p^^ ^ p^^) ^ 

where p = p^® p^ and Q*^^ (A^ > 0) denotes the A^-fold Kronecker power of Q. It is readily 
verified that the general formula for qubits is 

W(p) = Q^^0V'' P„,opPAf,{M-m) , (18) 

where p is any (not necessarily factorable) density matrix. Similarly, the inverse is given by 

U~\ct) = J]""^ P™,o(Q'''^0O■)PM,(M-„^) . (19) 



Evidently Q^^0W-i(cr) = 2-^W(Q'''^0 cr). 



III. REALITY CHECK 

The following properties of the real density matrix cr are worth noting explicitly: 
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• In addition to being real, it is nonsymmetric with one fixed element o"oo = 1; 

• It contains all the same information that the Hermitian density matrix does (since 
they are related by the bijection U). 

• It is diagonal if and only if the Hermitian density matrix is diagonal (which is why we 
defined an to be the coefficient of the diagonal Pauli matrix cr^). 

• It has the same tensor product structure as the Hermitian density matrix, since (as 
shown by Eq. (17)) U maps Kronecker products to Kronecker products. 

In addition to these nice analytic features, the real density matrix can also be quite useful for 
displaying the results of quantum state tomography: the determination of density matrices 
from experimental data. In most cases to date, the real or imaginary parts of the Hermitian 
density matrix have been displayed using two-dimensional bar graphs (see e.g (Nielsen and 
Chuang, 2000, §7.7.4)). Although useful, such a plot must both omit information and 
exhibit redundant information. Real density matrices are definitely superior in this respect, 
and sometimes may also exhibit the underlying symmetry of a state more clearly. Bar graphs 
of the real density matrix are shown below for both the diagonal Hilbert space basis as well 
as the Bell basis, illustrating how easily these states may be distinguished. Further examples 
with experimental NMR data may be found in Havel et al. (2002). 

The fact that the mapping between Pauli matrix coefficients and the entries of the Her- 
mitian density operator preserves the tensor product structure has been noted earlier by 
Pittenger and Rubin (2000) (and without doubt by many other researchers as well). In our 
present notation, their observation was based upon the following simple relation: 



O"00 CTlO 




1 1 




Poo Pio 






1 -1 




Pll Poi 



(20) 



While this is certainly a simpler relation than our operator sum, the "density matrix" on 
the right-hand side is not the usual Hermitian one, and the mapping between the two can 
be written explicitly only by using operator sums, supermatrices, or the like. Our goal here 
is to translate the usual operations and relations on Hermitian density matrices into the 
real domain, and the reordering of the entries of the real density matrix as above offers no 
advantage for this purpose. 
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FIG. 1 Plots of the diagonal Hilbert space basis (below) and of the Bell basis (above). All axes 
are dimensionless, and the labels on the horizontal axes correspond to the indices of the two-qubit 
real density matrix entries used in the main text. 



As our first example of such a translation, let us show how the usual criterion for the 
purity of the Hermitian density matrix can be carried over to the real domain: 

= ((Q 0Cr)> PA^,M-mPPm,0 ) 

\ 'm=0 / 

= (5^3^'^Pm,opPM,M-,n(Q®'^0O-)^) (21) 
= ((Q^"0^)(Q^"0^)^) 

= ((Q®^0Q^''0^)V) = (o-V) . 

In going to the last line, we have used the general relation ((A B) C^) = ((A C)^ B) for 
arbitrary conformant matrices A, B, C (Liitkepohl, 1996). This derivation easily generalizes 
to a formula for the ensemble-average expectation values of any observable with Hermitian 
matrix n and corresponding real matrix u = h({iJ,), showing that 

</.|p> = 2^(a.P> = (i^-'ct) = 2-^(^|a>. (22) 
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We close this section by noting that the partial trace operation corresponds simply to 
extracting a principal submatrix of the real density matrix (Somaroo et al., 1998). 



IV. LIFE IN THE REAL WORLD 

While the expectation values of observables carry over to the real domain without signif- 
icant complication, things become distinctly more challenging when it comes to integrating 
the equations of motion. In the case of a single qubit, it is readily verified that the commu- 
tator with an arbitrary Hamiltonian fi = U~^{v) becomes 

[<T, i/| = U[[p, ^i]) = zPoi(i^En<T - <TEni/)Poi 

(Til Z/io — CTlO J^ll 

(23) 

(Joi 1^11 - (711 Z/Ql O-ioZ/Ql - O-Qli^lO 

wherein the matrix entries are the components of the usual vector cross product. 

This equation of motion is most simply integrated by considering the matrix representa- 
tion of the commutation superoperator defined by which we henceforth assume without 
loss of generality has i^oo = 0. Letting | X) denote the column vector of height (M + 1)^ 
obtained by stacking the columns of the (M + 1) x (M + 1) matrix X on top of one another 
in left-to-right order, and applying the well-known identity 



I AXB ) = (B^ O A) I X ) 
(see e.g. (Liitkepohl, 1996)), we find that^ 



(24) 
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i/^Eii 



00 



(25) 



Since TZ^ = —WulflZi^ (II^P = 2(1/^1/) = this one-sided matrix differential 

equation is easily integrated by using the Cayley-Hamilton theorem to exponentiate the 



^ The factor of 1/2 here does not mean that the rotation is spinorial, but rather that the rate of rotation 
is I/V2 times the Hilbcrt-Schmidt norm of the Pauli matrix that generates it, while we pick up another 
factor of I/V2 on transforming to the real domain. 
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(lower-right 3x3 block of the) coefficient matrix (Najfeld and Havel, 1995), obtaining 

(26) 



K^)> - ( Poo ^ Poo + 7^. + ^ - ) nl ) I .(0) ) 



r(0) + zsin(||i.||t/2)[^(0), t>| - (l - cos(||i.|| t/2)) [[^(0), i> | , i> | ), 
wherein i> = and it is readily shown that 

= - ||i^inPoo®Poo - Eoo®Eoo) (27) 

is — lli^ll^ times the projection onto the plane orthogonal to the unit vector \u) = (^'1'''. 
The whole formula can thus be expressed more geometrically as 

\d-{t)) = (i>|d-) |i>) + sin(||i/||t/2) |i>) X Id") + cos(||i/||t/2) |i>) X (|i>) X (28) 

where & = cr(0) — Eqo and "x" is the cross product of the (last three components of the) 
vectors it connects. This is a standard expression for rotation of the three-dimensional vector 
I & ) about the axis 1 1>) by an angle ||i/||t/2. 

The extension of these formulae to general multi-particle commutators is not straight- 
forward, since the tensor product of commutation superoperators is not simply related to 
the tensor products of their underlying operators. Nevertheless, we can give a reasonably 
simple formula for the two-particle commutator with a factorable Hamiltonian, which is the 
most important case in practice. This is based upon a geometric algebra expression for the 
commutator of a tensor product of two three-dimensional vectors (or equivalently in the 
present context, traceless 2x2 Hermitian matrices), which is derived in Havel and Dor an 
(2002b): 

2[A®C,B®D] = (A|B)(Poo® [C,D]) + ([A, B] ® Pqo) ( C | D ) . (29) 

Letting a = W(A), etc. be the corresponding real matrices, this translates to: 

2[a®c, b®d| = (a|b)(Eoo® [c,d]) + ( [a, b] ® Eoo) ( c | d ) (30) 

This formula is easily extended to the case in which aoo, • • • , c^oo 7^ by multilinearity; in 
the following, however, we will need only the case in which aoo = cqo = 1, which introduces 
two additional terms: 

[Eoo®c, b®d| ^ i[Poo®C,B®D] = iB®[C,D] ^ib®|c,d]; ^^^^ 
[a®Eoo, b®d| ^ i[A®Poo,B®D] = i [ A, B] ® D ^ i |a, b] ® d . 
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Finally, we shall need the general result (Havel, 2003; Liitkepohl, 1996): 

|X®Y) = (Poo®/C22®Poo) (|X)® |Y)), (32) 
for 2x2 matrices X & Y, where the two-particle commutation matrix JC22 is given by 



3 3 



10 
10 
10 
1 



(33) 



We are interested in the case that a = cr^, b = i/^, c = cr^ and d = i/^, i.e. we have 
a factorizable two-particle state <t^ ® cr^ evolving under a bi-axial interaction ® u"^ . To 
express this more compactly, we define the matrix .Sj^i via 



' 




^01 
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^01 


















with an analogous definition in ( i/^ | <t^ ) | Eqo ) + | ) 
(31) & (32) give us 



I ) 



(34) 



5.2 I 0-2). Then Eqs. (25), (30), 



> 



= (Poo®/C22®Poo)|((^M'^')|Eoo>®|K, + (35) 

= (Poo ® /C22 ® Poo) i (5.1 ® 7?..2 + n^i ® 5.2) ( I ) (g) I cr^ )) 
= (Poo ® /C22 ® Poo) i (5.1 ® 7?..2 + 7?..i ® 5.2) ■ ■ ■ 

• • • (Poo ® /C22 ® Poo) I ® 0-2 ) . 

Since left or right multiplication of a 4D column or row vector by IZ^ gives the cross 
product of the last three components of that vector with [xio, Xoi, Xn], it may be seen that 
«S.i, IZjyi are mutually annihilating (i.e. «S.i7?-.i = 7?..i «S.i = 0), and similarly for «S.2, 
7?..2. As a result, the two terms on the last line of Eq. (35) commute and their exponential 
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factorizes. It is moreover easily shown that = ||x|| and TZ^. = — ||x|| 7?.^, so the 
overall integral is 

|cr(t)) = (Poo ® /C22 ®Poo) Exp(5^i ®7?.^2 t/4) Exp(7?.^i ® 5^2 t/4) ■•■ 

(36) 

■■■ (Poo ® /C22 ® Poo) |o-(0)), 
where cr(0) = cr^ a'^ and (letting Pqq be the 2x2 identity tensored with itself 4 times) 

Exp ( 5,1 ® n^2 t/A) = P|f + — ^ ' J ' ' (5,1 ® n^2) + 

11^ llll^ II 

l-cos(||i/i||i/2||t/4) ^ ^2 
11^ II 11^ II 

with an almost identical expression for Exp(7?.,i ® 5,2 t/A). The square in the last term 
may be evaluated by combining Eq. (27) with 

Si = |x)(x| + ||xf Eoo®Eoo . (38) 

Because 5,i®7?.,2 and 7?.,i®5,2 are mutually annihilating, the product of their exponen- 
tials expands to only five terms, two pairs of which have identical trigonometric coefficients. 
On pulling the right-hand column operator "| )" back to the left in Eq. (36), essentially 
reversing what we did to derive the differential version in Eq. (35), we obtain the integrated 
equation of motion we have been seeking: 



cr{t) = cr^ ® cr^ + sin(||i/^||||i/^|| t/4) [cr^ ® cr^, (g) 

(39) 

■■■ {l-cos{\\u'\\\y\\t/A))llc7'®a',u'®iy'' 



.2 ,.l^,.2T\ ^1 ^2 



Equations (30-31) tell us (more or less) what the geometric interpretation of the two particle 
commutator is, so we turn our attention to the double commutator. Geometric algebra shows 
that the double commutator of tensor products of the corresponding traceless 2x2 Hermitian 
matrices reduces to 

2 [[A® C, B^D], B^D] 

= [(A|B)(Poo® [C,D]) + ([A,B]®Poo)(C|D), B®D] (40) 
= (A|B)B® [[C,D],D] + [[A,B],B]®D(C|D) , 

so we will be done once we figure out what the double commutator of 2 x 2 matrices is. On 
expanding the commutator and using the fact that for such matrices the anticommutator 
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satisfies AB + BA = ( A | B ) Pqo, we get (including the real analogs): 

[[[a, bl,b]] ^ [[A, B],B] = AB2-2BAB + B2A 

= A||Bf - 2(BA + AB- AB)B (41) 
= 2(A||Bf - (A|B)B) ^ a||bf - (a|b)b . 

Back-substitution of this and the corresponding expression for [[C, D], D] into the preceding 
equation now yields: 

[[A(g)C, B® D], B®D] 

= (A|B)B(g) (C||Df - (C|D)D) + (A||Bf - (A|B)B) ®D(C|D) (42) 
= (A|B) ||Df (B® C) + (C|D) ||Bf (A®D) - 2 ( A | B ) ( C | D ) (B ® D). 

Since no matrix products occur in this expression, we may transliterate directly to the 
corresponding real expression including the additional terms arising from setting Ooo = 
Coo = 1 (cf. Eq. (31)): 

4[[(a+Eoo)®(c + Eoo), b®d|, b®d| = 

(a|b) ||df (b®c) + (c|d) ||b||2(a®d) ■■■ (43) 
- 2(a|b) (c|d) (b®d) ■■■ 
+ 2 [b® [[c, d] + [[a, bl ®d, b®d|. 

The first of these additional terms (given on the last line of the above equation) may be 
evaluated via Eqs. (30), (31) & (41) as indicated below, 

2[b®[[c,dl,b®d| = ||bf (Eoo®[[[c,dI,dI) + ([[b,bI®Eoo)(Ic,dl|d> 

(44) 

= ||bf (Eoo®(||df c - (c|d)d)) , 

with an analogous expression for the remaining term. We now put in our previous values 
a = <T^ = cr^ — Eoo, b = i/^, c = = cr^ — Eoo and d = and expand the result fully to 
get: 
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+ ||i/2||2(Eoo®o-2 + d-i ®Eoo) ■■■ (45) 

- (\\u'\\'{cT'\u')Eoo®iy' + {cT'\i.')\\u'fu'(^Eoo). 

Finally, we divide through by to get the normalized "vectors" and u"^, replace 

the cr's by a s inside the traces (which doesn't change their values) and recombine terms to 
obtain: 

{a' - {a'\i.')i.')0{{cr'\i.')i.' + Eoo) ■■■ (46) 

This has a fairly simple interpretation: The rejection of &^ from is tensored with the 
projection of onto (plus the usual scalar part of 1) and added to the projection of &^ 
onto (plus the scalar part) tensored with the rejection of ar^ from u^. This should be 
contrasted with the single commutator (Eq. (30)), wherein the inner and outer products of 
each pair are tensored together both ways and added. 



V. META-METAMORPHOSIS 

In the previous section we have shown how the Hamiltonians usually assumed for quantum 
computing with qubits can be integrated entirely within the real domain. For a single qubit, 
the results could be interpreted as a simple Bloch vector rotation. With a bi-axial interaction 
between two qubits, we also found that that the integrated expression had a reasonably 
nice geometric interpretation. Algebraically, however, it is usually easier to integrate in 
the Hermitian domain, simply because Hermitian matrices are easier to diagonalize. In 
this section, therefore, we shall derive formulae by which matrix representations of general 
superoperators can be translated from the Hermitian into the real domain, along with some 
specific examples of their utility. 

Using the identity given in Eq. (24) together with our operator sum expressions for U and 
its inverse (Eqs. (18) & (19)), it is straightforward to show that an arbitrary superoperator 
S with matrix representation <S acting on p via p i-^ «S | p ) transforms into the real domain 
according to 

. M \ / \ — 

S Q ( ^ Pm,(M-™') ® Pm',0 I^IY. PM,iM-m) ® Pm,0 ) Q , (47) 
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where Q = Diag(| Q®^)). The superoperator S may be written in operator sum form 
versus the basis of elementary matrices (Havel, 2003) as 

M,M M,M , M,M M,M 

^\P) = EE 4Efc.pE,A = E E 4(E£,®Efc,)|p). (48) 

i,j=0 k,e=o ' i,j=0 k,e=o 

On substituting the second of these equations (sans | p )) into the first and rearranging things 
a bit, the transformed superoperator becomes 

,M,M M,M M,M . 

2"'^S(EE4 E {'PM,iM-rn')'EijPMXM-rn))®{Prn'fi'EkiPrn,o)jQ. (49) 
^ i,j=0 k/=0 m,m'=0 ^ 

Unfortunately, because each factor in the above Kronecker product depends on both indices 
m and m', this cannot be regarded as a transformation of the elementary matrix basis into 
the real domain. 

Better insight can be obtained by looking at how the Choi matrix Choi{S) of the su- 
peroperator transforms. This may be obtained from the propagating matrix «S simply by 
replacing Kronecker products of the elementary matrices by dyadic products of the corre- 
sponding columnized basis, but in the opposite order (Havel, 2003). This task is facilitated 
by expressing the left- and right-multiplication by diagonal matrices as a Hadamard product, 
using the well-known formula (Havel et al., 2001) 

Diag(a)XDiag^(b) = (ab^) X , (50) 

which is essentially a special case of Eq. (32). This allows the transformed superoperator to 
be rewritten as 

Af,M M,M 

2-^(|Q^^>(Q^^|)0 5^ E 4 ■•■ 

i,j=0 k,e=o 
M,M 

■■■ E] ((PM,(M-m')^ejPM,{M-m)) ® {Pm'fi^kiPrnfifj ■ (51) 

m,m'=0 

The advantage of this form is that the Hadamard product commutes with the Choi operator 
(since it rearranges the entries of the product's operands identically), giving us 

M,M M,M 

Choi{S) ^ 2-^C/ioi(|Q®^>(Q^^|)0 5] E 4 

i,i=0 k,l=0 

M,M 

••■ E ^^<^*(^(PA^,(A^-m') Efj PA/,{A/-m)) ® (Pm',oEfcjPm,o)j 
m,m'=0 
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M,M M,M 
i,j=0 k,i=0 

M,M 

■■■ (\Pm',o'EikiPm,o){PM,[A'I-m')'EiejPM,{M-m)\J 
m,m'=0 

M,M M,M M,M 
jj'=0 fc,£=0 m,m'=0 



■ • ■ ( I (PM,(Af-m) ® PM,(Af-m')) (52) 

M,M MM M,M 

m,m'=0 i,j=0 k,£=0 



(Ejj Efcf) (PAf,(Af-m) ® PA/,(A/-m')) 
M' 

" ""Af',(Af'-m) ) 



2-^ (q""^® Q«^) E Pm.o Cho*(5) P 



m=0 

\2 1 _ r)2Ar 



where M' = {M + if - 1 = 2^^ - 1 and we have used the relation | Efc^ ) ( Eej | = (e^ 
ek){ej e^)"^ = E^j E^^. 

In other words, the Choi matrix of a superoperator maps into the real domain much like 
a density matrix on twice as many qubits. In fact we can write the real transformation 
matrix T, which acts on the real density matrix as T \(t), in the following compact form: 

r = 2-^ Choi(^U{Choi{S))^Q\ (53) 

Turning this around, we also find that we can express the Choi matrix of S in the Pauli 
basis as 

Choi{S) = U^\Choi{rQ^)) 

M,M M,M 

^ E E ^~'(4C/.oi((E,,0E,OQ^)) 



i,j=0 k,i=0 
M,M M,M 



E E [choi{E,j E,,,) Choi{ I 1 1^ )( Q®^0 Q^^ I ) 

i,j=0 k/=0 
M,M M,M 

= E E 4^-HE,,0E,,)0((Q«^0Q«^)0(ll^))) 

i,j=0 k,i=0 

M,M M,M 

= 2-^^((Q«^0Q«^) (11^)) E E 4 (p..® p..) , 

i,j=0 k,e=o 

where 1 denotes a column vector of ones of the appropriate size. 



(54) 
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It is time for our examples! We shall begin with operator sums for single qubit rotations, 
and go on to show how rotations about the axis as well as the crlcr'^ interaction between 
two qubits can be compactly described in the real domain using Hadamard products. We 
close by showing that this description also extends quite nicely to cr^ dephasing as well as 
nonunital relaxation back towards a nonrandom equilibrium state (i.e. T2 and Ti relaxation 
in NMR parlance). 

Consider the Choi matrix of the propagator which rotates a single qubit (Eqs. (25-27) 
above): 



CAioi(Poo®Poo + sin(||i/||t/2)7?.i> + (1 - cos(||i/|| t/2)) 72.^) 



1 








1 






























-z>n 


ho 










+ sin(||i/||t/2) 





















hi 


-ho 


1 








1 _ 




_ 


-1^01 1^01 



















- 





+ 



(55) 



;i-cos(||i.||t/2)) 



11 



hohi 

hohi 
'^01 - ^11 hohl hohi 



l^Ol'^U 



I'Oll^ll 



^01 - ^ 



10. 



where the "hat" on the i/'s indicates normalization by A general operator sum can be 
derived from this matrix (as well as by expanding the implied commutators via Eq. (23)), 
but since this is a bit involved we shall restrict ourselves to rotations by an angle d = 111^11^/2 
about the x, y or z coordinate axes. On substituting z>io = 1 and Pqi = 1^11 = into the 
above Choi matrix we obtain the following singular value decomposition: 

cos(??/2) sin(?9/2) 

sin(^9/2) -cos(i9/2) 
-sin(?9/2) cos(??/2) 

cos(i9/2) sin(^9/2) 

cos(i?/2) 10 

sin(i9/2) 1 



1 








cos('i9) 











sin('i9) 











— sin(^?) 


1 








cos('/9) 



(56) 



as may be readily verified using the usual half-angle formulae. The corresponding operator 
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sum for rotation by 1!} about the x-axis is simply: 



*W2)Pio 



cos(^9/2) 



sin(^9/2) 



cos(^?/2) -sin(i9/2) 
sin(^?/2) cos(?9/2) 

sin(^?/2) cos(^?/2) 
-cos(^?/2) sin(^9/2) 



1 CTOl 

CTlO 0-11 

1 O'Ol 

0"io o-ii 



+ ■•• 

1 
-1 



(57) 



In a similar fashion, it can be shown that the operator sum for a y-rotation is: 



«(i?/2)Poi 



cos(i9/2) 
sin(^9/2) 



1 


O'Ol 






0-11 




1 










-1 





cos(^?/2) -sin(?9/2) 
sin(^9/2) cos(^?/2) 

1 aoi sin(?9/2) cos(?9/2) 

aio ail -cos(^?/2) sm{^/2) 



(58) 



For a z-rotation, on the other hand, the Choi matrix turns out to be rank 4 with singular 
value decomposition: 



1 





cos('i9) 




1 








1 





— sin(^?) 










1 


-1 





sin(i9) 













1 


1 





cos(^9) 





1 




_ 1 








-1 




'cos2(^?/2) 























cos(^?/2) sin(i?/2) 

















cos(?9/2) sin(^9/2) 

sm\d/2) 

This corresponds to the operator sum 



1 








1 





1 


-1 








1 


1 





1 








-1 



(59) 



«W2)Pi 



(60) 



= cos2(i9/2)cr + sin2(i9/2)PncrPu + ■■■ 

•■■ zcos(?9/2)sin(^?/2)(Pioo-Poi + Poi^-Pio) , 
which has the pleasant feature that the trigonometric functions occur as scalar factors in 
each term and not embedded in the operators. This enables us to use Eq. (50) to rewrite it 
in terms of the Hadamard product as follows: 



*W2)Pnpg »W2)Pii\ 




1 cos('i9) 





1 


O"01 








cos('i9) 1 




0"lO 


o-ii 




1 






- sin(i9) 




1 








( 













1 




sin(i?) 






0-11 



+ 



1 

1 



(61) 
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As our next example, wherein the Hadamard product enables even greater simplifications, 
consider an "Ising-type" interaction between two qubits of the form cr^ (E) cr^ = P33, which 
is also known as "weak scalar coupling" in NMR (Ernst et al., 1987). An operator sum 
expression for this could be obtained by expanding the general formula given in Eq. (39), 
but because this Hamiltonian is again diagonal in the cr^ eigenbasis a simpler expression can 
be obtained directly starting from the diagonal matrix of the corresponding propagator, i.e. 



j-jP337rJt/2^giP337rJt/2 



,iP337rJt/2 



® e 



-lP33Tji/2 



)|P) 



Diag(|J(t)))|p) = Jit)\p) 



(62) 



where 



3{t) 



1 

-tirJt 



g mjt g iTT-Jt 



1 

-mjt 



Tjt 



-injt 



(63) 



1 1 

e 1 1 e 

2 g nrjt g nrjt 

On transforming this into the real domain via Eq. (47), we find that ^ has been converted 
into the sum of the diagonal matrix 

/ 



C{t) = Diag 



1 cos(7rJt) cos(7rJt) 1 
cos(7rJt) 1 1 cos(7rJt) 

cos(7rJt) 1 1 cos(7rJt) 



1 



cos(7rJt) cos(7rJt) 



(64) 



and the anti-diagonal matrix 

/ 



Pi5,o5(t) = Diag 



V 





sin(7rJt) 
sin(7r Jt) 




sin(7rJt) sin(7rJt) 





- sin(7r Jt) — sin(7r Jt) 





sin(7rJt) 
sin(7r Jt) 




(65) 



where the (self-inverse) left factor of Pi5,o 



simply reverses the order of the rows. The 



nonzero entries of the Choi matrix of C{t) +«S(t) turn out to comprise two 4x4 blocks along 
the diagonal, which are exactly the two matrices above, i.e. 



C{t) = Bisig{\£cChoi{C{t))£l)) 
and Pi5,o5(t) = maig{\esChoi{S{t))Sl)) 



(66) 
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wherein 



i=0 



67 ( 



and 



= ^ei(ei (8)e3_i) 



(67) 



i=0 



project out the rows / columns of their respective blocks. 

Thus we can obtain the desired operator sum representation by computing the eigenvalues 
and eigenvectors of the 4x4 symmetric matrices 



C(t) = ScChoi{C{t))Sl and S(t) = SsChoi{S{t))Sl 



(68) 



letting the operators' matrices be the diagonal / anti-diagonal matrices formed from the 
entries of these eigenvectors, and multiplying each term in the sum by the corresponding 
eigenvalue. The results are 



C(t)|cr) = 1(1 + cos(7r Jt)) cr + 1(1 — cos(7rJt)) P33 cr P 



33 



(69) 



and 



S{t)\a) 



sm{iTJt) P3oDiag([ 1, 1, 1, -1]) Diag([ 1, 1, 1,-1]) P 



30 



- P3oDiag([-l, 1, 1, l])o-Diag([-l, 1, 1, 1])P 



30 



(70) 



By using Eq. (50) to replace these operator sums by Hadamard products and taking advan- 
tage of the symmetry of S(t), however, we can obtain an even simpler expression, namely 



{C{t)+S{t))\(T) = |C(t)0O-- S(t)0(P3oO-P3o)> . 



(71) 



Finally, we show how one can also use Hadamard products with the real density matrix 
to describe simple relaxation processes, in a manner similar to that described in Havel et al. 
(2001) for the usual Hermitian density matrix. For a single qubit undergoing Ti (dissipation) 
and T2 (decoherence) relaxation, the time derivative is given by: 



dta{t) 



R (T{t) 



1/T2 





1 aoi(t) 


I/T2 1/Ti 




o'io(^) <yii{t) 



(72) 



Assuming that these relaxation processes are uncorrelated, this can immediately be extended 
to any number of qubits using the fact that the Hadamard product satisfies the mixed 
product formula with the Kronecker product (Eq. (15)). In the case of two qubits relaxing 
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with Hadamard relaxation matrices R^, R^, for example, we obtain 



(Ri®(ll^) 


1/T| 



® R2) a{t) 
1/T| 1/n 



1/Ti + 1/T| 
1/Tf i/n + l/r| l/T^ + l/Tf 
l/T^i + 1/T| 1/Ti 1/T,i + 1/T| 

i/Tl + i/r| i/T^i + i/Tf i/Tl + i/r| i/r^i + i/rf 



o-(t) . 



(73) 



where 1 is a 4 x 1 vector of I's. The fact that uncorrelated Ti as well as T2 relaxation can 
be extended so easily to multiple spins in this way is actually a significant advantage of 
the real density matrix over the Hermitian, since in the latter case the diagonal terms are 
mixtures of terms decaying at differing rates, substantially complicating their treatment via 
Hadamard products (Havel et al., 2001). 

When correlations are present, however, these advantages are largely lost, since then the 
off-diagonal entries of the real density matrix consist of mixtures of terms with differing 
decay rates (fortunately, Ti relaxation is usually largely uncorrelated (Ernst et al., 1987)). 
Let us work through the case of two qubits in detail, assuming for simplicity that the T2 
relaxation processes at the two qubits are totally correlated and have the same rate I/T2, 
as for example in an NMR gradient-diffusion experiment (Havel et al., 2001). In this case 
the Hadamard relaxation matrix for the Hermitian density matrix has the form 



R = — 



1 

f2 



114 
10 1 
10 1 
4 110 



(74) 



and the corresponding 16 x 16 diagonal relaxation superoperator 



7^ = Diag(|R)) 



(75) 



is easily exponentiated into a diagonal matrix of survival probabilities for the entries of the 
(traceless part of the) Hermitian density matrix. In this case, however, it turns out to be 
almost as easy, but more revealing, to convert 7?. into the real domain and perform the 
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integration there. The result of the first step is 
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i > 
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¥2 
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.0 












































0. 



(76) 



which again has nonzero entries only on its diagonal and its anti-diagonal. Representing the 
action of this matrix on | cr) as a sum of a Hadamard product and a Hadamard product 
coupled with row/column inversion as in Eq. (71), we obtain the real equation of motion: 



-T2dt(T{t) 



112 
10 2 1 
12 1 
2 110 



o-(t) + 



-2 

2 

2 

-2 



(P30O-(t)P30) (77) 



It may readily be verified that the operations on cr which occur in the two terms of this 
expression commute, and hence this equation can be integrated by exponentiating them 
separately. With the first term this leads to a simple Hadamard (entrywise) exponential 
(Havel et al., 2001), namely 



/ 



D(t) = ExPq 



t 

7^ 



v 



1 1 
1 
1 
2 



2 
2 



2 
1 
1 



1 


e-t/T2 


e-t/T2 


g-2t/72 




1 


g-2t/r2 






g-2t/T2 


1 




g-2i/T2 


e-t/T2 


e-t/T2 


1 



(78) 



1 1 J / 

To exponentiate the second term, we note that the Hadamard product is with 2 P03 
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2 cr2 (g) cr2 and resort briefly to superoperators in order to simplify the exponential as follows: 
Exp( -(2t/T2) Diagd Po3 )) (P30 ® P30) ) 
= Er=o ( " ^2^/^^) I )) (P30 ® P30) ) 7 A:! 



Poo® Poo + Diag(|P3o))^ 



(2t/T2 



fc=l 



2k\ 



Diagd Po3))(P 



30 



30. 



E 

fc=0 



(2t/T2 



\2fe+l 



(79) 



(2A;+ 1)! 



00 



00 



Diagd P3o))(cosli(2t/T2) "1) ■■■ 

- Diag( I Po3 )) (P30 ® P30) sinh(2t/r2) . 

This derivation relies on the facts that P30 ® P30 squares to the identity, Diag(|Po3)) 
squares to Diag(| P30)), and each commutes with the other. 

Going back to operator sum notation and abbreviating cr = cr(0), we thus obtain in all: 

cr(t) = D(t) © (o- + (cosh(2t/T2) - 1) P30 cr ■ ■ ■ 

- sinh(2t/T2) Po3 © (P3ocrP3o)) 



O"00 

cToi e" 



(Toi e 
0-11 



0-026 



-t/T2 



^03, (T30) 



A(t; 0-12, -0-21) 0-136" 



~t/T2 



-t/T2 



A{t; 0-21, -0-12) 



0-22 



A{t; 0-30, 0-03) 0-31 e" 



0-326 



-t/T2 



0-236 
0-33 



(80) 



where A(t; x, y) = (cosh(2t/T2)x + sinh(2t/T2) y)) exp(— 2t/T2). From this we see that 
the anti-diagonal entries decoher into mixtures with their symmetrically placed opposites 
in the real density matrix. These mixtures correspond to the real and imaginary parts of 
the P12 = P21 entries in the Hermitian density matrix, otherwise known as zero-quantum 
coherences, which are immune to correlated noise (Ernst et al., 1987). 



VI. EPILOGUE 

We have seen that one can, with some effort, do pretty much everything with the real 
density matrix that one could with the usual Hermitian one. This may be useful as a didactic 
device, or in calculations with experimental (e.g. NMR) data where it is desirable to keep 
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the experimentally measured values of the observables in sight at all times. This work is 
also a good demonstration of the power of Choi matrix decompositions as a means of finding 
operator sum representations of linear superoperators (Havel, 2003). 

Although the Hermitian density matrix is expected to be better suited, by and large, for 
the purposes of numerical calculations, it is worth emphasizing that for theoretical and/or 
expository purposes the compact but lucid notation of geometric algebra offers significant 
advantages over any matrix formalism. In this regard, we point out that Havel et al. (2003) 
have recently introduced a parity- even (rather than reverse-even, aka Hermitian) multi- 
qubit density operator via geometric algebra, which generalizes the multi-particle space-time 
algebra introduced for isolated systems to open multi-qubit systems. It is our hope that 
in due course such a geometric formulation may provide new insights into some of the 
conceptual problems that underlie quantum physics. 

The existence of the real density matrix is further of some theoretical interest, since 
it provides a coordinate ring within which one can study the issues of entanglement and 
decoherence via invariant theoretic methods (Grassl et al., 1998; Makhlin, 2002). There are 
intimate connections between invariant theory and geometric algebra, and it is often easier 
to automate symbolic computations in an invariant ring than it is at the more abstract level 
of geometric algebra (Havel, 1997, 2001; Sturmfels, 1993). 
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